function [X_dot] = plantAll(x, Fd, Md, F_plant, M_plant, q_i_d, q_i, a_all, ui_n)

global g m_o Jx Jy Jz m Quadrotor_Num
global a1 a2 a3 l_i
global b1 b2 b3
global kt kd 
X_dot = zeros(14,1);
X_plant = zeros(12, Quadrotor_Num);
Ux = zeros(4,1);
Uy = zeros(4,1);

for i = 2 : Quadrotor_Num + 1
    Ux(i) = cos(x(1,i))*sin(x(3,i))*cos(x(5,i))+sin(x(1,i))*sin(x(5,i));
    Uy(i) = cos(x(1,i))*sin(x(3,i))*sin(x(5,i))-sin(x(1,i))*cos(x(5,i));
end


e3 = [0; 0; 1];

% 这里的x应该是货物的x，暂且先这样书写，后面要改 TODO
X_dot(1) = x(2,1);
X_dot(2) = Md(1) / Jx + a1 * x(4,1) * x(6,1);
X_dot(3) = x(4,1);
X_dot(4) = Md(2) / Jy + a2 * x(2,1) * x(6,1);
X_dot(5) = x(6,1);
X_dot(6) = Md(3) / Jz + a3 * x(2,1) * x(4,1);
X_dot(7) = x(8,1);
X_dot(8) = 1/m_o * Fd(1);
X_dot(9) = x(10,1);
X_dot(10) = 1/m_o * Fd(2);
X_dot(11) = x(12,1);
X_dot(12) = 1/m_o * Fd(3) + g;

% 绳索的角速度
for plant_i = 1 : Quadrotor_Num
    X_s_dot(1,plant_i) = x()
end
X_dot(13) = x(14);
X_dot(14) = 1/l_i * skew(q_i) * a_all - 1/m/l_i * skew(q_i) * ui_n;

% 无人机的位姿信息
for plant_i = 2 : Quadrotor_Num + 1
    X_plant(1,plant_i) = x(2,plant_i);
    X_plant(2,plant_i) = x(4,plant_i) * x(6,plant_i) * a1 + M_plant(1,plant_i) * b1;
    X_plant(3,plant_i) = x(4);
    X_plant(4,plant_i) = x(2,plant_i) * x(6,plant_i) * a2 + M_plant(2,plant_i) * b2;
    X_plant(5,plant_i) = x(6,plant_i);
    X_plant(6,plant_i) = x(2,plant_i) * x(4,plant_i) * a3 + M_plant(3,plant_i) * b3;
    
    % 位置状态
    X_plant(7,plant_i)  = x(8,plant_i); 
    X_plant(8,plant_i)  = Ux(plant_i) * F_plant(plant_i)/m;
    X_plant(9,plant_i)  = x(10,plant_i); 
    X_plant(10,plant_i) = Uy(plant_i) * F_plant(plant_i)/m;
    X_plant(11,plant_i) = x(12,plant_i); 
    X_plant(12,plant_i) = -g + cos(x(1,plant_i)) * cos(x(3, plant_i)) * F_plant(plant_i)/m;
end


end 


function S = skew(v)
    % v 为 3x1 向量, S 则为 3x3 反对称矩阵,满足 S*w = v x w
    S = [   0,   -v(3),  v(2);
         v(3),     0,   -v(1);
        -v(2),   v(1),    0   ];
end